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HIGHLIGHTS 


•  Integrated  model:  a  phase-field  model  and  an  electrochemical  cell  one  are  linked. 

•  Using  the  integrated  model,  we  study  the  role  of  Ni  coarsening  on  SOFC  performance. 

•  Simple  scaling  laws  describe  the  evolution  of  key  microstructural  features. 

•  Ni  particle  size  strongly  affects  morphological  evolution  in  Ni— YSZ  SOFC  anodes. 

•  Integrated  approach  facilitates  identification  of  parameters  associated  with  stable  microstructures  and/or  optimal  performance. 
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In  solid  oxide  fuel  cells  (SOFCs),  Ni  coarsening  in  porous  anodes  that  are  comprised  of  Ni  and  yttria 
stabilized  zirconia  (YSZ)  leads  to  changes  in  several  microstructural  attributes,  which  affect  the  elec¬ 
trochemical  performance.  Herein  we  present  an  integrated  modeling  approach,  where  a  dynamic 
mesoscale  phase  field  model  is  linked  with  a  stationary  macroscale  electrochemical  cell  level  model  in 
order  to  assess  the  role  of  Ni  coarsening  on  the  performance  of  SOFCs.  The  phase  field  model  is  capable  of 
capturing  the  morphological  evolution  of  Ni  and  accounting  for  its  polycrystalline  nature,  while  the 
electrochemical  model  encompasses  the  entire  set  of  processes  of  gas  transport,  electronic  and  ionic 
conduction  as  well  as  the  electrochemical  reactions.  Microstructural  features  are  extracted  from  the 
phase  field  model  as  anode  systems  evolve  over  time  and  employed  as  effective  properties  in  the 
electrochemical  model.  Simulation  results  highlight  the  importance  of  Ni  and  YSZ  particle  size  and  ratio 
on  both  the  microstructural  stability  and  electrochemical  performance  of  SOFCs.  In  particular,  it  is  shown 
that,  for  the  classes  of  microstructures  employed  in  this  work,  coarsening  of  Ni  particles  can  either 
improve  or  diminish  the  maximum  power  density  relative  to  the  as-sintered  ones,  depending  on  the 
initial  particle  size. 

©  2013  Elsevier  B.V.  All  rights  reserved. 


*  Corresponding  author.  Department  of  Mechanical  and  Aerospace  Engineering, 
Princeton  University,  Princeton,  NJ  08544,  USA.  Tel.:  +1  609  2589126;  fax:  +1  609 
2585877. 

E-mail  addresses:  fabdelja@princeton.edu  (F.  Abdeljawad),  benvoelker@gmail. 
com  (B.  Volker),  rsdavis@princeton.edu  (R.  Davis),  rmcm@engineering.ucsb.edu 
(R.M.  McMeeking),  mhaataja@princeton.edu  (M.  Haataja). 

0378-7753 /$  -  see  front  matter  ©  2013  Elsevier  B.V.  All  rights  reserved. 
http://dx.doi.Org/10.1016/j.jpowsour.2013.10.121 


1.  Introduction 

Solid  oxide  fuel  cells  (SOFCs)  are  regarded  as  alternative  devices 
for  cleaner  and  more  efficient  energy  conversion  [1—3].  With  the 
advantages  of  high  efficiency  of  60%-80%  [4,5]  and  reduced  carbon 
emissions  [6],  SOFCs  are  being  developed  for  various  stationary  and 
auxiliary  commercial  power  applications  with  a  desired  lifetime  of 
40,000  h  or  more  [7].  In  conventional  SOFCs,  which  operate  at 
temperatures  in  the  range  800-950  °C  [8,9],  a  thin  film  of  yttria- 
stabilized  zirconia  (YSZ)  is  utilized  as  the  oxygen  ion  conducting 
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ceramic  membrane.  A  porous  composite  of  lanthanum  strontium 
manganate  (LSM)  and  YSZ  is  typically  used  for  the  cathode  and  a 
state-of-the-art  anode  is  a  porous  Ni-YSZ  composite  10]. 

In  SOFCs,  oxygen  is  electrochemically  reduced  at  the  cathode, 
from  which  the  oxygen  ions  subsequently  diffuse  through  the 
electrolyte  and  react  with  the  oxidized  fuel  at  the  anode.  The  anode 
in  turn  is  a  porous,  complex  three  phase  material,  within  which  all 
phases  (pore,  conducting  metal,  and  ceramic  oxide)  percolate.  The 
contiguity  of  the  pore  phase  ensures  effective  transport  of  the  fuel 
to  the  catalytic  sites,  while  contiguous  metallic  and  ceramic  oxide 
phases  serve  as  conduits  of  the  charges  (electrons  and  oxygen  ions) 
involved  in  the  electrochemical  reactions  taking  place  in  the  vi¬ 
cinity  of  three-phase  boundaries  (TPBs),  which  constitute  the 
electrochemically  active  regions  [11].  The  anode  must  maintain  a 
stable  microstructure  over  extended  periods  of  time  at  elevated 
temperatures,  as  well  as  high  electron  and  ion  conductivity  and 
catalytic  activity. 

Electrochemical  performance  degradation  in  SOFC  cermet  an¬ 
odes  is  usually  attributed  to  several  phenomena  and  processes  at 
the  microstructural  level,  such  as  poisoning  by  fuel  impurities, 
coking,  thermal  and  redox  cycling  [11-16  ,  as  well  as  microstruc¬ 
tural  coarsening  processes  associated  with  Ni  particles  in  Ni-based 
cermet  anodes  [9,17,18].  Coarsening  is  a  competitive  growth  pro¬ 
cess,  which  leads  to  a  reduction  in  the  total  interfacial  energy. 
Spatial  variations  in  the  mean  interface  curvature  induce  diffu- 
sional  fluxes,  which  transfer  material  from  regions  of  high  mean 
curvature  to  ones  associated  with  lower  mean  curvature.  In  the  case 
of  isotropic  interfacial  energies,  coarsening  leads  to  a  reduction  in 
the  total  interfacial  area  [19,20]. 

Indeed,  Ni  coarsening  in  SOFC  anodes  has  been  the  subject  of 
several  experimental  investigations  in  recent  years  [11,21—24]. 
These  studies  have  revealed  that  performance  degradation  can  be 
directly  linked  to  changes  in  several  microstructural  features 
associated  with  coarsening  processes,  such  as  extent  of  TPBs,  total 
interfacial  areas,  and  contiguity  of  anode  phases.  In  particular,  the 
density  of  TPB  lines,  which  is  often  employed  as  a  simple  metric  to 
predict  electrochemical  performance,  is  found  to  decrease  due  to  Ni 
coarsening  [23].  Simwonis  et  al.  [21]  studied  Ni-YSZ  composite 
anodes  under  operating  environments  and  observed  a  decrease  in 
the  electrical  conductivity  due  to  Ni  coarsening.  Moreover,  Ni  phase 
coarsening  is  found  to  affect  phase  contiguity  [25].  A  recent  study 
by  Cronin  et  al.  [26]  revealed  that  Ni  coarsening  affects  several 
microstructural  attributes,  such  as  pore  percolation,  and  total  pore 
and  Ni-pore  interfacial  areas. 

While  the  aforementioned  experimental  studies  have  provided 
important  insights  to  the  evolution  of  anode  microstructures  in 
operating  environments,  they  do  not  examine  the  effects  of  inter¬ 
related  physical  phenomena  on  coarsening  rates.  Computational 
models  on  the  other  hand  are  ideally  suited  for  exploring  the 
various  hypothesized  mechanisms  that  control  microstructural 
evolution  processes  in  SOFC  anodes.  Existing  computational  studies 
of  Ni  coarsening  in  SOFC  anodes  include  two-particle  systems  27], 
phase  field  models  [28,29],  and  multi-scale  frameworks  [30].  Our 
approach  builds  on  and  significantly  extends  these  previous 
computational  studies  by  incorporating  the  coupling  between 
microstructural  evolution  and  electrochemical  performance  while 
accounting  for  the  polycrystalline  nature  of  the  Ni  phase. 

In  this  work,  we  focus,  in  a  parametric  study,  on  the  role  of 
anode  microstructure,  characterized  by  Ni  particle  size  and  Ni  to 
YSZ  particle  size  ratio,  on  the  topological  evolution  of  SOFC  anodes 
and  its  impact  on  electrochemical  performance.  More  specifically, 
we  present  an  integrated  modeling  approach  to  examine  the 
impact  of  Ni  phase  coarsening  in  Ni-YSZ  porous  anodes  on  the 
performance  of  SOFC  cells.  In  this  approach,  several  representative 
anode  microstructures  with  tunable  features  (such  as  phase 


volume  fractions  and  particle  sizes)  are  first  generated  “virtually”. 
Next,  a  mesoscale  diffuse-interface  model,  which  is  capable  of 
capturing  Ni-YSZ  interaction,  and  accounting  for  the  poly¬ 
crystalline  nature  of  Ni,  is  employed  to  quantitatively  investigate 
the  coarsening  of  Ni  phase  particles  in  said  microstructures.  Then, 
morphological  features  are  extracted  from  the  phase  field  model  as 
the  anode  systems  evolve  over  time  and  employed  as  effective 
input  parameters  to  a  recently  developed  macroscale  cell  level 
model  to  quantitatively  evaluate  the  time-dependent  electro¬ 
chemical  performance  [31  ].  Here,  we  note  that  while  the  electro¬ 
chemical  model  itself  is  time-independent,  its  input  parameters, 
which  represent  effective  properties  that  depend  on  the  anode 
microstructure,  are  evolving  in  time  and  tracked  via  the  phase  field 
model. 

Starting  from  a  binary  mixture  of  numerically  sintered  Ni  and 
YSZ  particles,  our  simulation  results  conclusively  demonstrate  that 
the  overall  TPB  length  decreases  monotonically  during  coarsening. 
Furthermore,  the  reduction  in  the  overall  TPB  length  obeys  a  simple 
scaling  law  in  time.  In  addition,  the  contiguity  of  the  Ni  phase, 
which  also  monotonically  decreases  during  coarsening,  displays  an 
interesting  size  effect,  such  that  systems  with  initially  small  Ni 
particles  lose  contiguity  faster.  Finally,  the  overall  electrochemical 
performance  can  display  a  non-monotonic  time-dependence  dur¬ 
ing  coarsening,  and  even  improve  over  time.  For  initially  coarse 
microstructures,  the  electrochemical  performance  is  dominated  by 
the  electrochemical  reactions  taking  place  in  the  vicinity  of  TPBs, 
the  length  of  which  decreases  monotonically.  For  initially  fine  mi¬ 
crostructures,  on  the  other  hand,  electrochemical  performance  is 
dominated  by  gaseous  transport  processes  through  the  pore  phase, 
leading  to  an  increase  in  the  maximum  power  density  over  time  as 
the  average  pore  size  increases  due  to  coarsening. 

More  broadly  speaking,  simulation  results  highlight  the 
importance  of  initial  anode  microstructure  on  the  stability  and 
performance  of  SOFC  anodes.  Our  modeling  approach  can  be  used 
as  a  design  tool  to  identify  combinations  of  morphological  pa¬ 
rameters,  such  as  Ni  and  YSZ  particle  size  and  ratio  that  yield 
porous  anode  systems  with  optimal  microstructural  stability  and 
electrochemical  performance  characteristics. 

The  rest  of  the  manuscript  is  organized  as  follows.  In  Section  2 
we  first  outline  both  the  continuum  modeling  framework  and  the 
cell  level  approach  to  evaluate  the  electrochemical  performance 
based  on  microstructural  parameters.  Then,  in  Section  3  we 
describe  how  the  initial  microstructures  are  generated  and  char¬ 
acterized.  Quantitative  results  from  coarsening  studies  and  their 
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Fig.  1.  A  block  diagram  illustrating  the  integrated  modeling  approach  that  is  used  to 
examine  the  role  of  Ni  coarsening  in  Ni-YSZ  porous  anodes  on  the  electrochemical 
performance  of  SOFCs.  The  arrow  outlined  in  blue  represents  the  feedback  from  the 
electrochemical  model  to  the  input  microstructures,  which  in  turn  closes  the  design 
iteration  loop.  This  step  is  not  implemented  in  the  work  presented  here.  (For  inter¬ 
pretation  of  the  references  to  color  in  this  figure  legend,  the  reader  is  referred  to  the 
web  version  of  this  article.) 
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implications  to  electrochemical  performance  are  in  turn  discussed 
in  Section  4.  Finally,  concluding  remarks  are  presented  in  Section  5. 


2.  An  integrated  modeling  approach 

Our  computational  approach  links  a  mesoscale  phase  field 
model  for  Ni  coarsening  in  Ni-YSZ  porous  anodes  to  a  macroscale 
cell  level  model,  which  evaluates  the  electrochemical  performance 
of  SOFCs.  The  aforementioned  models  are  implemented  as  follows. 
First,  the  phase  field  model  is  used  to  track  the  evolution  of  several 
microstructural  features  as  “virtual”  anode  systems  evolve  over 
time.  Next,  these  features  are  extracted  from  the  phase  field  model 
and  used  as  inputs  to  the  cell  level  model  in  order  to  evaluate  the 
electrochemical  performance  of  these  systems.  A  block  diagram 
illustrating  our  modeling  approach  is  shown  in  Fig.  1. 

Below,  we  first  describe  the  salient  features  of  the  phase  field 
formalism  employed  in  our  coarsening  studies.  Then,  an  overview 
of  the  cell  level  approach  to  predict  the  electrochemical  perfor¬ 
mance,  for  a  given  microstructure,  is  provided.  Finally,  the  nu¬ 
merical  generation  and  microstructure  characterization  of  “virtual” 
anode  systems  are  described. 


2.1.  The  phase  field  formalism 


The  challenge  in  microstructural  evolution  studies  resides  in 
accurately  capturing  the  complex  topological  changes  (such  as 
merging  of  two  or  more  particles  or  the  disappearance  of  one)  that 
take  place  during  the  coarsening  process.  The  so-called  phase  field 
approach,  which  has  already  been  employed  in  studies  of  SOFC 
microstructures  by  other  groups  [28-30],  elegantly  incorporates 
such  topological  changes  in  a  numerically  convenient  manner  as 
discussed  next. 

The  starting  point  of  our  phase-field  model  is  the  introduction  of 
order  parameters  for  the  phases  in  porous  Ni— YSZ  anodes.  Addi¬ 
tional  structural  order  parameters  are  also  introduced  to  represent 
the  polycrystalline  nature  of  the  Ni  phase  and  account  for  interfa¬ 
cial  energies  associated  with  Ni  grain  boundaries  (GBs).  Next,  a  total 
free  energy  T tot  that  describes  the  state  of  the  anode  system  is 
constructed  in  terms  of  the  order  parameters  and  the  dynamics  of 
the  system  follow  from  the  minimization  of  Consistent  with 
experimental  observations  [22,26],  a  static,  or  non-evolving,  order 
parameter  0(x)  is  assigned  to  the  YSZ  phase,  reflecting  its  micro- 
structural  stability.  In  contrast,  the  evolution  of  the  Ni  phase  is 
tracked  with  the  aid  of  an  order  parameter  c(x,t)  in  combination 
with  structural  order  parameters  {^/(x,t),  i  =  1,...,  n^}  that  describe 
the  shape  and  crystallographic  orientation  of  the  Ni  grains.  Flere,  nv 
is  the  number  of  structural  order  parameters  needed  to  resolve  the 
texture  of  the  Ni  phase.  The  order  parameters  are  chosen  such  that 
(a)  in  the  YSZ  phase  ( d,c )  =  (1,-1 ),  (b)  in  the  Ni  phase  (0,c)  =  (0,  +1 ), 
and  (c)  in  the  pore  space  (0,c)  =  (0,-1).  Furthermore,  0(x)  e  (0,1) 
and  c(x,t)  e  (-1,1 )  define  YSZ  and  Ni  interfaces,  respectively.  For  the 
structural  order  parameters,  the  equilibrium  values  are  chosen  such 
that  within  non-Ni  phases  (YSZ  and  pore):  {77*  =  0,  i  =  1,...,  nv }  and 
within  a  Ni  grain  {ip  =  1,  =  0,  ij  =  1,. As  is  shown  in  Fig.  2, 

order  parameters  take  on  a  constant  value  within  the  bulk  phases 
and  vary  rapidly-but-smoothly  across  interfaces. 

As  in  standard  Ginzburg-Landau  formalism  and  motivated  by 
the  work  of  Wang  [32]  on  solid-state  sintering,  the  following  total 
chemical  free  energy  is  employed  for  an  SOFC  anode  system: 


F tot[c,  0,  {Vi}] 
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Fig.  2.  A  schematic  illustrating  the  order  parameters  used  in  the  phase  field  model.  6 
and  c  parameters  are  for  the  YSZ  and  Ni  phases,  respectively.  r)m  and  r?n  are  structural 
order  parameters  for  Ni  grains  m  and  n,  respectively. 


where  w  and  Iq  are  so-called  gradient  energy  coefficients  employed 
to  tune  the  Ni  surface  and  grain  boundary  energies.  The  first  term  on 
the  right  hand  side  (/bulk)  is  a  fourth-order  Landau  polynomial  that 
defines  the  homogeneous  bulk  free  energy  associated  with  the  Ni 
phase  and  multiple  Ni  grains.  In  particular,  /bulk  is  constructed  such 
that  the  equilibrium  values  for  the  order  parameters,  i.e.  9/buik/ 
6c  =  0/buik/<h7z  =  0,  are  (c  =  -1,  rp  =  0,  i  =  1,...,  nv)  within  non-Ni 
phases  and  (c  =  +1,  rn  =  +1,  ijj^i  =  0,  ij  =  1,...,  nv)  within  a  Ni 
grain.  A  convenient  form  for/bUik  that  satisfies  these  requirements  is 

/buik(c,H-})  =  mH(c2  -  l)2  +Ah|T(c+  l)2  +  3(1 

~2(3-c)f>?  +  3 

i  =  1 

(2) 

where  mu  and  AH  are  parameters  we  employ  to  tune  the  interfacial 
energies  and  interface  widths.  In  this  work,  it  is  assumed  that  GB 
interfaces  are  isotropic,  yielding  Iq  =  k.  The  last  term  (fwe tt)  on  the 
right  hand  side  of  Eq.  (1)  accounts  for  the  interaction  between  the 
Ni  and  YSZ  phases,  and  the  resulting  wetting  angle  at  triple  junc¬ 
tions.  In  this  work,  /wett  is  written  as 

/wett(G^)  —  —  ?lc|v0|  +  —  C^/3^  ,  (3) 

where  and  £2  are  model  parameters  employed  to  tune  the 
equilibrium  Ni— YSZ  wetting  angle,  experimentally  observed  to  be 

>117°  [33]. 

With  the  proposed  total  chemical  free  energy  functional  in  Eq. 
(1)  and  following  the  treatment  of  Cahn  and  Hilliard  [34],  the 
interfacial  energies  and  widths  can  be  evaluated  analytically.  In 
particular,  we  obtain  the  Ni  surface  energy  ys  ^4w^2mH/3 
+ky/AH/3  and  interfacial  width  <5s~2w/^2mH,  while  Ni  GB  energy 
7gb  -2/<^Ah/3  and  its  width  <5cB~k/v^H-  Therefore,  the  phase 
field  model  parameters  w,  mn,  Ah  and  k  are  uniquely  determined 
from  the  Ni  phase  interfacial  and  GB  energies  and  widths. 

Within  the  Ginzburg-Landau  formalism  and  with  the  aid  of 
variational  derivatives,  the  governing  equations  for  the  evolution  of 
c(x,t)  and  {^/(x,t),  i  =  nv}  for  (x,t)  e  V  x  [0,T]  c  u3xU ,  where  V 

is  the  domain  volume  and  T  is  the  final  time,  can  be  stated 
respectively  as 
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2.2.  A  cell  level  model  for  the  electrochemical  performance  of  SOFC 
(4a)  systems 


at 


tot 

'rsp 


i  =  1 ,2,-,^, 


(4b) 


where  M  denotes  the  Ni  phase  mobility  function  and  Li  is  a  kinetic 
parameter  associated  with  Ni  GB  mobility,  where  under  the  isot¬ 
ropy  assumption  of  GB  interfaces,  L;  =  L.  At  SOFC  operating  tem¬ 
peratures,  Ni  bulk  self  diffusion  ~10-16  m2  s-1  [35],  whereas  Ni 
surface  self  diffusion  DNj  ~  1CT11  m2  s'1  [36].  We  also  note  that  the 
Zr4+  surface  diffusion  coefficient  in  YSZ  is  ~10^18  m2  s_1  [37], 
justifying  the  assumption  that  YSZ  particles  do  not  significantly 
evolve  during  Ni  coarsening.  Since  Ni  diffusivity  along  YSZ  in¬ 
terfaces  is  not  well  characterized,  it  is  assumed  that  the  coarsening 
of  Ni  in  Ni-YSZ  cermet  anodes  is  dominated  by  Ni  surface  diffusion 
along  Ni-pore  interfaces.  In  this  work,  Ni  phase  mobility  has  the 
general  form  M  =  MNiQ(c,Vc,#,V<9),  where  Mm  is  the  Ni  phase  atomic 
mobility,  and  the  detailed  form  of  the  dimensionless  function  Q  is 
given  below  in  Eq.  (7). 

The  spatio-temporal  evolution  equations  for  c(x,t)  and  r/j(x,t)  are 
made  non-dimensional  by  introducing  non-dimensional  spatial 
coordinates  x  =  x/d0  and  time  t  =  t/i,  where  d0  and  t  denote  the 
characteristic  length  (voxel  edge  length)  and  time  scales,  respec¬ 
tively.  For  convenience,  we  also  introduce  the  non-dimensional 
quantities  w2  =  mH/(2E),  mH  =  mHd2/(<52E),  =  AHdl/(5^E), 

?i  =  ?!/(E<52),  and  f2  =  ?2^o/(^s)»  where  E  is  a  reference  energy 
density.  Ni  mobility  can  be  written  as  MNi  =  DNi/(dp/dc)Ceq  = 

Dni/  (92/buik/^^2)ceq*  Therefore,  the  characteristic  time  t  can  be 
expressed  as 


Due  to  the  coarsening  of  the  Ni  phase,  the  microstructure  in  the 
reactive  anode  layer  evolves  with  time.  In  order  to  systematically 
capture  the  effects  of  these  microstructural  changes  on  the  elec¬ 
trochemical  performance  of  the  SOFC  system,  we  employ  a  nu¬ 
merical  cell  level  model  encompassing  the  entire  set  of  processes  of 
gas  transport,  electronic  and  ionic  conduction  as  well  as  the  elec¬ 
trochemical  reactions.  For  completeness,  a  brief  overview  of  the 
approach  is  provided  next;  the  reader  is  referred  to  Ref.  [31  for  a 
more  detailed  exposition  of  the  framework. 

The  model  geometry  consists  of  a  porous  anode  support  and 
reactive  layer,  a  dense  electrolyte  layer,  and  porous  cathode  reactive 
and  current  collector  layers  as  illustrated  in  Fig.  3.  In  this  work,  we 
only  consider  the  effects  of  the  temporal  microstructural  changes 
associated  with  Ni  coarsening  in  the  reactive  anode  interlayer, 
while  the  effective  properties  of  all  other  layers,  i.e.,  anode  support, 
dense  electrolyte  and  cathode  remain  constant  over  time. 

In  our  approach,  we  extract  representative  3D  anode  micro¬ 
structures  from  the  phase-field  model  at  regular  time  intervals,  and 
compute  the  effective  microstructural  parameters  that  characterize 
the  electrochemical  performance  of  the  reactive  anode  layer.  These 
effective  properties  are  then  transferred  into  the  stationary  one¬ 
dimensional  cell  level  model  and  serve  as  input  parameters  for 
computing  the  electrochemical  performance.  This  approach  yields 
temporal  snapshots  of  the  electrochemical  performance,  which 
enables  us  to  quantitatively  examine  the  role  Ni  phase  coarsening 
on  the  evolution  of  electrochemical  performance  of  SOFCs. 

Now,  the  equations  for  mass  transport  and  migration  of  elec¬ 
trons  and  oxygen  ions  in  the  reactive  anode  layer  take  the  form  of 
Poisson  equations: 


mH^(92/bulk/9C2)^  16W2d4 
£Dn  i^s 

For  the  systems  we  examine  in  this  work,  where  d0  ~  60  nm, 
5S  ~  1  nm,  and  DN i  ~  1CT11  m2  s-1,  the  resulting  characteristic  time 
t  ~  40  s.  As  for  the  governing  equations  for  the  non-conserved 
parameters,  we  introduce  the  non-dimensional  L  =  xLEd^/d2 

~2  i 

and  k  =  k2/(EdQB).  GB  mobility  Mqb  and  energy  ygb  are  related  to 
the  model  parameters  through 


Lk2 


MgbYgb 


Ife2d2^Ik2DNi62 
T  16w2d2 


(6) 


where  Eq.  (5)  has  been  employed  in  the  last  step. 

The  resulting  non-dimensional  spatio-temporal  evolution 
equations  for  c(x,t)  and  rn(x,t )  have  the  same  form  as  the  original 
ones.  For  the  remainder  of  this  paper,  we  drop  the  tildes  on  all 
parameters  for  notational  convenience.  Finally,  the  form  that  is 
used  for  the  Ni  phase  mobility 


M  =  MNi  ^1.0  -  tanh^|V(9|2//i1^tanh^|Vc|2//i2^) 

+  M*tanh(|V0|2/h1)]],  (7) 

where  c  =  tanh (c/h3).  The  model  parameters  hi,  h2  and  h 3  are 
employed  to  tune  the  range  of  c  and  6  values  that  contribute  to  the 
mobility  function,  while  M*  defines  a  mobility  contribution  for  Ni 
diffusing  along  YSZ-pore  interfaces.  In  this  work,  it  is  assumed  that 
diffusion  along  Ni-pore  interfaces  is  the  dominant  mass  transport 
mechanism. 


VNql  (Rpore?  ^00 ) 


Areacdan 


2  F 


(8) 


—  —  Areacd'an,  (9) 

v(-4sfZ^YSz)  —  —  Areact7arb  (10) 

where  Na  denotes  the  rate  of  mass  transport  for  gas  species  a  and  is 
a  function  of  the  average  pore  radius,  Rp0re,  and  the  pore  tortuosity 
factor,  a Further,  Areact  is  the  reactive  TPB  length,  ian  the  local 
Faradaic  current  density,  and  F  the  Faradaic  constant.  The  diffusive 
gas  transport  is  modeled  using  the  dusty-gas  approach  [38  ,  which 
accounts  for  both  molecular  and  Knudsen  diffusion,  and  detailed 
expressions  used  for  Na(Rpo re,a°o)  can  be  found  in  Ref.  [31].  Ohm’s 
law  describes  the  transport  of  electrons  and  ions,  where  (7^f  is  the 
effective  electronic  conductivity  of  the  Ni  phase,  the  effective 
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Fig.  3.  A  schematic  illustration  of  an  anode  supported  SOFC  cell  that  is  used  for  the 
analysis  in  this  work. 
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Table  1 

Parameters  used  in  the  ID  cell  model  for  the  functional  layers:  anode  support  (AS),  reactive  anode  interlayer  (AN),  reactive  cathode  interlayer  (CA)  and  cathode  current 
collector  (CC).  Note  that  effective  properties  of  the  reactive  anode  interlayer  that  evolve  due  to  Ni  coarsening  are  designated  “evolving”  in  this  table. 


Parameter 

AS 

AN 

CA 

CC 

Ref. 

L  [pm] 

1000 

20 

20 

50 

[39,56] 

Po [Pa] 

101,300 

101,300 

101,300 

101,300 

[17] 

(0Ni»  0YSZ) 

(0.26,  0.26) 

(0.30,  0.45) 

(0.38,  0.38) 

(0.55,  0.0) 

[39] 

0pore 

0.48 

0.25 

0.24 

0.45 

[39] 

K ysz  [pm] 

1.0 

0.5 

1.0 

(l?Ni»  ^pore)[pm] 

(1.0,  1.62) 

Evolving 

(0.5,  0.29) 

(1.0,  1.45) 

[39,56] 

Areact  [mm”3] 

Evolving 

5.72  x  1012 

CCoo 

1.65 

Evolving 

2.2 

1.71 

(Teff 

Ni 

0.033 

Evolving 

0.059 

0.214 

^YSZ 

0.059 

0.059 

06  0r 

0.5 

0.5 

0.5 

0.5 

[57] 

MHz  [mol  1  kg] 

2.016  x  10-3 

2.016  x  10-3 

MH2o  [mol1  kg] 

M0z  [mol”1  kg] 

18.02  x  10-3 

18.02  x  10-3 

32.00  x  10-3 

32.00  x  10-3 

MN2  [mol1  kg] 
(£h2.£h2o)[I<] 

(59.7,  809.1) 

(59.7,  809.1) 

28.01  x  10-3 

28.01  x  10-3 

[58,59] 

(£o2.  £n2) [K]  o 

(106.7,71.4) 

(106.7,  71.4) 

[58,59] 

(Fh2,  FH2o)  [A] 

(2.827,  2.641) 

(2.827,  2.641) 

[58,59] 

(f02,  rNz)  [A] 

(3.467,  3.798) 

(3.467,  3.798) 

[58,59] 

Eh2  LI  mol  '] 

120  x  103 

120  x  103 

[60] 

Eq2  [J  mol1] 

130  x  103 

130  x  103 

[60] 

Gh2  [J  mol-1] 

-147.1  x  103 

-147.1  x  103 

[61] 

CH2o  [J  mol-1] 

-454.9  x  103 

-454.9  x  103 

[61] 

Cq2  [J  mol-1] 

-228.0  x  103 

-228.0  x  103 

[61] 

G02  [J  mol-1] 

-236.4  x  103 

-236.4  x  103 

-236.4  x  103 

-236.4  x  103 

[61] 

ionic  conductivity  of  the  YSZ  phase,  $Ni  the  electric  potential  in  the 
Ni  phase  and  $Ysz  the  ionic  potential  in  the  YSZ  phase.  The  charge 
transfer  rate  correlated  with  the  hydrogen  oxidation  reaction  tak¬ 
ing  place  in  the  reactive  anode  layer  is  described  by  a  Butler— 
Volmer  equation  as 


(11) 


where  ft  and  are  the  forward  and  reverse  charge  transfer  co¬ 
efficients,  R  is  the  universal  gas  constant,  T  is  the  temperature  and 
£an  is  the  anodic  overpotential.  For  a  detailed  expression  of  the 
anodic  reference  exchange  current  density,  ian(ran),  where  the 
parameter  ran  is  a  fitting  constant  in  our  model,  the  reader  is 
referred  to  Ref.  [31].  Likewise,  we  describe  all  physical  processes  of 
gas  transport,  electron  and  ion  migration  and  electrochemical 
redox  reactions  in  the  remaining  layers:  electron  transport  and 
diffusion  of  hydrogen  and  water  vapor  in  the  anode  support  layer, 
oxygen  ion  transport  through  the  dense  electrolyte  layer,  the  oxy¬ 
gen  reduction  reaction  and  ion  transport  in  the  reactive  cathode 
layer  as  well  as  the  transport  of  electrons  and  the  diffusion  of  ox¬ 
ygen  and  nitrogen  (atmospheric  air)  through  the  cathode  current 
collector  and  reactive  cathode  layer. 

Our  model  [31]  is  parameterized  based  on  experimental  results 
reported  by  Zhao  and  Virkar  [39],  who  consider  a  laboratory  anode- 
supported  SOFC  button  cell.  Their  cell  is  operating  at  a  constant 


temperature  of  Tre f  =  1073  I<  with  hydrogen  and  H2O  vapor  partial 
pressures  at  the  anode  support  inlet  of  pHl  =  0.97p0,an  and 
Ph2o  =  Po,an  -Ph2*  and  oxygen  and  nitrogen  partial  pressures  at 
the  cathode  current  collector  inlet  of  p0l  =  0.21p0,ca  and 
Pn2  =  Po,ca  -  Po2  (atmospheric  air)  at  a  uniform  total  gas  pressure 
ofpo.an  =  Po.ca  =  101.3  kPa.  As  schematically  illustrated  in  Fig.  3,  the 
cell  consists  of  an  anode  support  layer  of  dimensions  Las  =  1000  pm, 
a  reactive  anode  interlayer  (Lan  =  20  pm),  a  dense  electrolyte  layer 
(Liyte  =  8  pm),  a  reactive  cathode  interlayer  (Lca  =  20  pm)  and  a 
current  collector  (Lcc  =  50  pm).  These  dimensions  and  the  model 
parameters  of  all  layers  are  chosen  to  correspond  to  the  experi¬ 
mental  values  published  by  Zhao  and  Virkar  [39],  and  summarized 
in  Table  1.  The  bulk  conductivities  of  Ni,  YSZ  and  LSM  at 
Tref  =  1073  K  are  taken  from  the  literature  [40]  as 

0Ni,o  —  3.03  x  104  Q  1  m  \  oysz.o  =  2.26  Cl  ^  m  ^  and 
0lsm,o  =  1-28  x  104  Q_1  m-1,  and  the  fitting  constants  of  the  model 
are  determined  as  Tan  =  4  x  10~4  A  m-1  and 

rca  =  3.75  x  10-4  A  m-1  by  using  a  best-fitting  approach.  Assuming 
a  steady  state  operating  condition  for  the  fuel  cell,  i.e.,  treating  gas 
flow  as  well  as  electronic  and  ionic  currents  as  time-independent, 
the  model  is  numerically  solved  using  the  commercial  finite- 
element  software  COMSOL  Multiphysics  3.5a,  with  cell  level 
boundary  conditions  as  described  in  Ref.  [31]. 


3.  SOFC  electrodes 


Table  2 

Properties  of  virtual  SOFC  anode  systems  considered  in  this  study.  The  as-prepared 
anode  active  layer  is  idealized  as  a  binary  mixture  of  electron  and  ion  conducting 
spherical  particles  with  radii  RN j  and  RY sz.  respectively.  For  all  cases,  the  volume 
fractions  of  Ni  (0Ni),  YSZ  (</>Y sz)  and  pore  (0pore)  phases  are  30  ±  2%,  45  ±  2%,  25  ±  2%, 
respectively. 


System  name 

^Ni/^YSZ 

Number  of 

Ni  particles 

Number  of 

YSZ  particles 

R1 

1.00 

60 

90 

R1.5 

0.67 

206 

92 

R2 

0.50 

503 

95 

R2.5 

0.40 

1010 

97 

3.1.  Numerical  generation  of  SOFC  electrode  microstructures 

Virtual  random  packing  systems  of  spherical  particles  are  utilized 
in  the  integrated  modeling  approach  for  the  microstructural  evolu¬ 
tion  of  SOFC  anodes  and  its  impact  on  the  electrochemical  perfor¬ 
mance.  SOFC  electrodes  (anode  and  cathode)  are  idealized  as  binary 
mixtures  of  electron  and  ion  conducting  particles.  The  influence  of  the 
local  anode  microstructure  characterized  by  the  Ni  to  YSZ  particle  size 
ratio  on  the  coarsening  of  Ni  and  the  subsequent  electrochemical 
performance  is  investigated.  Flere  it  is  assumed  that  cathode  does  not 
undergo  any  evolution  processes  that  influence  its  microstructure. 
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Fig.  4.  Perspective  views  of  the  virtual,  as-prepared  anode  microstructures.  Ni  and  YSZ  particles  are  represented  in  orange  and  blue,  respectively,  (a)  R1  system:  RNi  =  RY sz  (b)  R1.5 
system:  RNi  =  0.67RYsz  (c)  R2  system:  RNi  =  0.5RYsz  (d)  R2.5  system:  RNi  =  0.4RYSZ.  In  all  systems,  the  volume  fractions  of  Ni  (0Ni),  YSZ  (0YSZ)  and  pore  (0pore)  are  30%,  45%,  25%, 
respectively,  with  ±2%  added  to  all  fractions.  Note  that  the  YSZ  particle  size  is  fixed  in  all  systems  and  only  the  number  and  size  of  Ni  particles  are  varied.  (For  interpretation  of  the 
references  to  color  in  this  figure  legend,  the  reader  is  referred  to  the  web  version  of  this  article.) 


We  use  the  random  close  packing  algorithm  of  Gan  et  al.  [41,42] 
to  generate  binary  mixtures  of  spherical  particles  with  periodic 
boundary  conditions.  Then,  in  order  to  numerically  mimic  the  effect 
of  sintering  and  enable  the  formation  of  connected  particle  clusters, 
overlap  between  the  particles  is  generated  by  expanding  the  radii 
of  all  spheres  uniformly  and  proportional  to  the  radius  while 
keeping  the  center  of  each  sphere  at  a  fixed  position,  until  a  desired 
densification  is  reached.  Four  virtual  Ni-YSZ  anodes  are  examined, 
denoted  Rl,  R1.5,  R2,  R2.5,  with  Ni  to  YSZ  particle  radius  ratio  R^/ 
Rysz  =  1.0,  0.67,  0.50  and  0.40,  respectively.  For  all  systems  studied 
in  this  work,  the  resulting  volume  fractions  of  Ni  (0nO,  YSZ  (0Ysz) 
and  pore  (0pore)  phases  are  30  ±  2%,  45  ±  2%,  25  ±  2%,  respectively, 
similar  to  experimentally  characterized  Ni-YSZ  based  cermet  an¬ 
odes  22,23].  Furthermore,  in  this  parametric  study,  the  size  of  YSZ 
particles  is  kept  constant.  For  a  given  SOFC  anode  system,  or  a  given 
Rm/Rysz  ratio,  the  size  and  number  of  Ni  particles  are  varied  in  order 
to  attain  the  target  volume  fraction.  Table  2  lists  the  Ni  to  YSZ 
particle  size  ratio  and  the  number  of  Ni  and  YSZ  particles  for  each 
SOFC  anode  system  used  in  this  work.  Perspective  views  of  the  Rl, 
R1.5,  R2  and  R2.5  systems  depicting  Ni  and  YSZ  particles  and  the 
resulting  virtual  anode  microstructures  are  in  turn  shown  in  Fig.  4, 
where  Ni  and  YSZ  particles  are  colored  in  orange  and  blue, 
respectively.  Finally,  to  account  for  sample-to-sample  variations, 
coarsening  simulations  on  four  ensembles  for  each  SOFC  anode 


system  are  performed  and  results  are  averaged  for  the  use  of  the 
electrochemical  performance  analysis. 

3.2.  Microstructural  characterization 

Several  microstructural  attributes,  whose  evolution  directly  af¬ 
fects  the  electrochemical  performance,  are  monitored  during  the 
course  of  Ni  phase  coarsening  in  the  reactive  anode  layer.  Of  special 
interest  are  the  effective  material  properties  used  in  Eqs.  (8-10)  of 
the  cell  level  model:  the  reactive  TPB  Areact>  the  average  pore  radius 
Rpore*  and  the  pore  tortuosity  factor  as  well  as  the  effective 
conductivities  a^f  and  The  fields  c(x,t)  and  0(x,t)  are  extracted 
from  the  phase  field  model  at  several  time  steps  during  the  tem¬ 
poral  evolution  of  anode  systems,  and  used  to  describe  the  topology 
of  three-phase  anodes. 

3.2.1.  Percolation  and  three-phase  boundary 

We  define  the  percolation  degree  Pd  as  the  fraction  of  a  phase 
that  belongs  to  spatially  contiguous  clusters,  i.e.  percolating  clus¬ 
ters  that  span  the  simulation  box  from  end  to  end.  For  example, 
Pd  =  0.5  indicates  that  50%  of  the  Ni  phase  belongs  to  Ni  clusters 
that  are  spatially  contiguous.  Spatially  isolated  clusters  are  identi¬ 
fied  using  the  labeling  algorithm  developed  by  Floshen  and 
Kopelman  [43].  It  is  worth  mentioning  here  that  the  YSZ  phase  in  all 
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anode  systems  used  in  this  study  is  a  non-evolving  one  with  a 
constant  percolation  degree  ~99.5%.  Therefore,  it  is  sufficient  to 
only  monitor  the  temporal  evolution  of  the  Ni  phase  percolation 
degree  Pd. 

Junctions  in  the  microstructure  of  the  porous  anode,  where  the 
Ni,  YSZ  and  pore  phases  meet  define  the  three-phase  boundary, 
where  electrochemical  reactions  take  place.  We  distinguish  be¬ 
tween  the  total  TPB,  Atot,  which  accounts  for  all  junctions  where 
the  three  phases  meet,  and  the  reactive  TPB,  Areact»  which  only 
accounts  for  the  fraction  of  junctions  where  the  constituents  of  all 
three  phases  belong  to  percolated  clusters.  Using  the  percolation 
information  we  obtain  during  the  evolution  of  SOFC  anodes,  we 
determine  the  length  of  the  reactive  TPB,  Areact*  by  counting  voxel 
edges  where  percolated  Ni,  YSZ  and  pore  voxels  meet.  The  length  of 
reactive  TPB,  Areact,  in  porous  SOFC  electrodes  represents  an 
important  parameter  for  the  electrochemical  performance  as  it 
dictates  the  overall  charge  transfer  rate  of  the  redox  reactions. 

3.2.2.  Pore  space  properties 

Serving  as  input  parameters  for  modeling  the  diffusive  gas 
transport  through  the  pore  network  (cf.  Eq.  (8)),  the  average  pore 
radius,  Rpore,  and  the  pore  tortuosity  factor,  aM,  are  extracted  from 
the  voxel  discretized  electrode  microstructures.  In  order  to  deter¬ 
mine  ftpore,  we  first  define  the  pore  phase  field  2(x,t),  where  2  =  1 
denotes  the  pore  phase  and  zero  elsewhere.  Then,  we  compute  the 
two-point  correlation  function  g2(r,t)  =  ([2(r,t)  -  (2(r,t))] 
[S(r')  -  (2(r,t))])/g,  where  /i  =  (E2(r,t))  -  (2(r ,t))2  is  the  variance. 
The  average  pore  radius  Rpore  is  defined  as  the  full  width  at  half 
maximum  ofg2(r,t)  [44]. 

When  considering  the  diffusion  of  fuel  gas  through  porous 
media,  the  irregular  geometry  of  the  pore  space  shows  a  significant 
effect:  gas  molecules  are  forced  to  follow  tortuous  pathways 
through  the  pore  space,  which  slows  down  the  rate  of  diffusional 
transport.  In  the  dusty-gas  model  used,  this  effect  is  taken  into 
account  by  the  pore  tortuosity  factor,  aM.  A  Lattice-Monte-Carlo 
(LMC)  method  is  used  to  obtain  from  the  voxel  discretized 
electrode  microstructures.  A  large  number  of  so-called  tracers 
(typically  106)  are  placed  on  pore  voxels  and  perform  random 
walks.  During  each  LMC  step,  the  mean  square  displacement  of  all 
tracers  is  monitored.  Sufficient  LMC  steps  are  reached  once  the 
temporal  derivative  of  the  mean  square  displacement  becomes 
independent  of  time,  what  we  typically  see  after  2  x  106  LMC  steps. 
Then,  the  tortuosity  factor  is  determined  by  comparing  the  mean 
square  displacement  that  is  calculated  in  the  pore  space  geometry 
to  the  mean  square  displacement  of  unbiased  tracers  in  free  space. 

3.2.3.  Electron  and  oxygen  ion  migration 

As  embodied  in  Eqs.  (9)  and  (10),  Ohm’s  law  describes  the 
transport  of  oxygen  ions  through  the  YSZ  phase  as  well  as  the 
transport  of  electrons  through  the  Ni  phase,  where  the  effective 
conductivities  and  are  required  as  input  parameters.  A  fully 
detailed  description  of  our  approach  for  computing  the  effective 
conductivities  of  the  Ni  and  YSZ  phases  from  the  voxel  discretized 
microstructures  can  be  found  in  Ref.  [45].  This  approach  is  again 
based  on  the  Lattice-Monte-Carlo  (LMC)  method  and  utilizes 
percolation  information  obtained  in  the  cluster  identification  study. 

4.  Results  and  discussion 

4.1.  Role  ofNi  coarsening  on  microstructural  features  of  SOFC 
anodes 

Parameters  of  the  phase  field  model  are  chosen  such  that  they 
mimic  the  behavior  of  real  SOFC  anodes.  First,  the  parameters  of  the 
mobility  function  are  chosen  as  Mm  =  1.0,  M*  =  0.01,  h\  =  0.01, 


h2  =  0.005  and  /13  =  0.20,  which  ensures  that  transport  mechanisms 
controlling  Ni  phase  coarsening  is  dominated  by  Ni  diffusion  along 
Ni-pore  interfaces.  Young’s  equation  relates  surface  and  grain 
boundary  energies  to  the  dihedral  angle  at  junctions  where  grain 
boundaries  meet  free  surfaces.  At  SOFC  operating  temperatures,  the 
ratio  of  Ni  grain  boundary  to  surface  energy  7gb/Ts  =  0.25-0.5  [46], 
which  yields  dihedral  angles  ranging  between  150°  and  165°. 
Therefore,  in  the  phase  field  model,  we  set:  k  =  3.0,  w  =  V2  and 
m\\  =  1.0  and  AH  =  1.0,  which  yields  tgb/Ts  ~  0.4.  Ni  phase  grain 
boundary  mobility  Mgb  is  temperature-dependent  and  can  vary  by 
orders  of  magnitude  depending  on  the  geometrical  properties  of 
the  grain  boundary  and  impurity  content  47,48].  In  this  work,  the 
kinetic  coefficient  L  =  1.0,  which  corresponds  to  grain  boundary 
mobility  MGb  ~  10-11  m4  s-1  J  1  typically  found  in  metals  [49,50],  is 
employed.  Finally,  the  model  parameters  describing  Ni-YSZ  con¬ 
tact  angle  were  set  to  ?  1  =  5.0  and  ?2  =  1.0. 

To  verify  that  the  model  incorporates  the  appropriate  interfacial 
energetics,  we  carried  out  simulations  of  neck  growth  of  two 
originally  spherical  Ni  particles  and  wetting  properties  of  a  Ni 
particle  deposited  on  YSZ.  The  results  are  shown  in  Fig.  5(a)  and  (b), 
respectively.  As  can  be  seen  in  Fig.  5(a),  the  contact  between  two  Ni 
particles  at  the  end  of  the  simulation  is  characterized  by  a  dihedral 
angle  — 155°.  Furthermore,  as  shown  in  Fig.  5(b),  the  observed  Ni- 
YSZ  wetting  angle  is  -123°.  Both  angles  are  representative  of 
experimentally  observed  ones  in  the  Ni-YSZ  SOFC  anode  system. 

Flaving  demonstrated  that  the  model  properly  incorporates 
interfacial  energetics,  we  start  our  exploration  of  the  model  by 
examining  the  virtual  SOFC  anode  systems  with  various  Ni  to  YSZ 
particle  size  ratios,  or  RuilRysz-  Each  virtual  SOFC  cermet  anode  is 
discretized  to  produce  a  simulation  box  of  size  120  x  120  x  120 
voxels,  where  the  voxel  size  is  Ax  =  Ay  =  Az  =  1.0.  Periodic 
boundary  conditions  in  all  three  dimensions  are  used  in  the  phase 
field  simulations.  First,  we  qualitatively  investigate  the  morpho¬ 
logical  evolution  due  to  Ni  coarsening.  For  the  R2.5  anode  system, 
where  RNi  =  0.4RYsz»  panels  (a)-(c)  of  Fig.  6  depict  the  evolution  of 
Ni  phase  spatially  isolated  clusters,  where  every  Ni  cluster  is 
assigned  a  unique  color  at  (dimensionless)  coarsening  times  0, 103 
and  104,  respectively.  Initially,  a  single  and  highly  intertwined 
cluster  constitutes  most  of  the  Ni  phase.  As  the  anode  system 
evolves  over  time,  parts  of  the  Ni  skeleton  break  into  smaller 
clusters,  as  can  be  seen  in  Fig.  6(b)  and  (c),  which  eventually  leads 
to  loss  of  percolation  degree  Pd.  This  loss  of  percolation  degree  by 
fragmentation  of  large  clusters  into  smaller  ones  will  be  discussed 
in  more  detail  later  in  this  section. 


Fig.  5.  (a)  Simulation  of  the  morphological  evolution  of  two  Ni  grains  (orange)  where 
Tgb/Ts  ~  0.4  and  the  resulting  dihedral  angle  of  ~  155°.  The  dashed  blue  line  is  drawn 
to  indicate  the  spatial  location  of  the  grain  boundary,  (b)  Simulation  of  Ni  (orange) 
interaction  with  a  rigid  YSZ  (gray)  substrate  and  the  resulting  wetting  angle  of  ~  123°. 
In  both  panels,  dashed  and  solid  black  lines  represent  initial  and  equilibrium  config¬ 
urations,  respectively.  (For  interpretation  of  the  references  to  color  in  this  figure 
legend,  the  reader  is  referred  to  the  web  version  of  this  article.) 
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Fig.  6.  For  the  R2.5  anode  system  (RNi  =  0.4FYsz)-  (a)— (c)  Iso-surfaces  of  spatially  isolated  Ni  clusters  with  each  cluster  assigned  a  unique  color  at  time  (a)  t  =  0  (b)  t  =  103,  and  (c) 
t  =  104.  (d)— (f)  Temporal  evolution  of  the  Ni  phase  grain  morphology,  with  Ni  grain  boundaries  colored  in  black,  (d)  t  =  0,  (e)  t  =  103,  and  (f)  t  =  104.  Note  that  the  static  YSZ  particles 
are  omitted  from  the  visualizations  for  clarity. 


Next,  we  examine  the  temporal  evolution  of  the  Ni  phase 
morphology.  Panels  (d)-(f)  of  Fig.  6  depicts  the  temporal  evolution 
of  the  grain  morphology  of  the  Ni  phase  for  the  R2.5  anode  system. 
Initially,  a  collection  of  small  spherical  particles  constitutes  the  Ni 
phase.  As  the  anode  system  evolves  over  time,  Ni  grain  growth  is 
observed,  where  the  average  size  of  grains  increases  and  the 
number  of  Ni  grains  shrinks.  In  addition  to  Ni  coarsening,  where  the 
total  interfacial  energy  is  minimized,  grain  growth  leads  to  reduc¬ 
tion  in  the  total  interfacial  energy  associated  with  grain  boundaries. 

We  now  turn  our  attention  to  several  microstructural  features 
affecting  the  electrochemical  performance.  We  first  start  by 
examining  the  role  of  microstructure  (Ni  to  YSZ  particle  size  ratio) 
on  the  evolution  of  TPB  regions.  Within  the  phase  field  framework, 
TPB  regions  are  defined  as  spatial  regions  where  interfaces  of  all 
phase  fields  meet  and  overlap.  Therefore,  the  total  number  of  vol¬ 
ume  elements  (voxels),  where  the  Ni,  YSZ,  and  pore  phases  meet,  is 
regarded  as  a  measure  of  total  TPB,  or  A^ot.  Fig.  7(a)  is  a  plot  of  Af0t 


Fig.  7.  Temporal  evolution  of  (a)  total  TPB  A(ot  (b)  reactive  TPB  A*eact.  In  both  figures, 
TPB  is  measured  by  the  voxel  count  of  elements  that  are  spatially  located  at  Ni-YSZ- 
pore  interfaces. 


for  all  systems.  It  can  be  clearly  seen  that  for  a  fixed  YSZ  particle 
size,  systems  with  smaller  initial  Ni  particles  lead  to  higher  A£ot.  A 
quantity  of  critical  importance  when  modeling  electrochemical 
reactions  is  reactive  TPB,  or  A*eact,  which  is  defined  in  this  work  as 
the  total  number  of  elements,  where  percolating  Ni  clusters  meet 
the  YSZ  and  pore  phases.  Fig.  7(b)  is  a  plot  of  the  temporal  evolution 
of  A*eact.  At  late  coarsening  times,  systems  with  smaller  initial  Ni 
particles  experience  a  sharper  decrease  in  A*eact.  This  is  due  to  the 
loss  of  percolation  of  some  Ni  clusters,  as  will  be  shown  later  in  this 
section. 

It  is  worth  noting  that  the  TPB  data  is  expressed  in  non- 
dimensional  units  in  order  to  allow  for  the  investigation  of  elec¬ 
trochemical  performance  of  porous  anodes  with  a  broad  range  of  Ni 
and  YSZ  particle  size  and  ratio.  The  results  shown  in  Fig.  7  can  be 
made  dimensional  by  assigning  a  value  to  the  characteristic  length 
scale  d0,  which  sets  Ni  and  YSZ  particle  size  along  with  the  char¬ 
acteristic  time  t.  For  example,  assuming  d0  =  40  nm,  which  corre¬ 
sponds  to  t  ~  35  s,  a  simulation  box  of  (4.8  pm)3,  and  YSZ  particle 
size  of  1  pm  (RY sz  =  0.5  pm),  Table  3  lists  reactive  TPB  density  for 
two  SOFC  anode  systems,  R1  and  R2,  at  the  initial  configuration  and 
end  of  simulation,  which  corresponds  to  a  (dimensional)  coars¬ 
ening  time  ~  100  h.  Simulation  results  for  both  absolute  value  and 
relative  change  of  reactive  TPB  density  in  these  systems,  Table  3,  are 
in  good  agreement  with  experimental  observations,  where  a  factor 
of  ~2  reduction  in  reactive  TPB  density  was  observed  [25,26]. 

Table  3 

Assuming  a  characteristic  length  scale  dQ  =  40  nm,  reactive  TPB  density  for  two  SOFC 
anode  systems,  R1  and  R2,  with  Ni  and  YSZ  particle  radii  (RNi,  RY sz)  =  (0-5  pm, 
0.5  pm)  and  (0.25  pm,  0.5  pm),  respectively. 


RNi  [pm] 

Rysz  [pm] 

TPB  [mm"3] 

Initial 

After  ~100  h 

Change [%] 

0.5 

0.5 

5.24  x  1012 

3.04  x  1012 

-42 

0.25 

0.5 

5.63  x  1013 

1.76  x  1013 

-69 
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Fig.  8.  (a)  On  a  iog-iog  scale,  total  TPB  A£ot  vs.  coarsening  time,  (b)  The  pre-factor 
/  =  /(^Ni/^Ysz).  where  is  assumed  to  follow:  Atot  =  f  (Rm/Rvsz)tp ■  The  dashed 
line  is  drawn  for  clarity. 


Now,  a  close  examination  of  Fig.  7(a)  reveals  a  self  similar 
behavior  for  the  temporal  evolution  of  AJot  during  Ni  coarsening. 
Fig.  8(a)  is  a  log-log  plot  of  AJot  vs.  time,  where  A£ot  is  shown  to  be 
consistent  with  a  power  law  of  the  form 

Atot  =fxtp,  (12) 

where  the  exponent  p  =  -0.16  ±  0.02  and  /  =  /{RniIRysz)  is  a  pre¬ 
factor,  which  encompasses  the  effects  of  anode  microstructure. 
Fig.  8(b)  in  turn  demonstrates  the  dependency  of  the  pre-factor/on 
the  Ni  to  YSZ  particle  size  ratio  RniIRysz,  where  the  dashed  line  is 
drawn  as  a  guide  to  the  eye. 

Next,  we  examine  the  evolution  of  pore  space  during  the 
microstructural  evolution  of  SOFC  anodes.  Fig.  9(a)  is  a  plot  of  the 
average  pore  radius  Rpore  normalized  by  the  YSZ  particle  radius 
Rysz,  which  is  constant  (YSZ  is  a  non-evolving  phase)  and  the  same 
in  all  systems.  In  all  SOFC  anode  systems,  Rpore  increases  during  the 
evolution  of  Ni  phase.  We  also  notice  that  SOFC  anode  systems 
having  larger  initial  Ni  particles  yield  larger  pores  Rpore-  Fig-  9(b)  is  a 
plot  of  the  total  pore  space  interfacial  area  (with  both  Ni  and  YSZ) 
Apore  normalized  by  its  initial  value  Apore  (t  =  0).  All  SOFC  anodes 
experience  a  continuous  decrease  in  pore  interfacial  area  with 
systems  having  smaller  initial  Ni  particles  experiencing  a  larger 
drop  in  pore  interfacial  area.  At  the  end  of  coarsening  simulations, 
the  total  decrease  in  pore  interfacial  area  -25%,  40%,  50%  and  55% 
for  R1  (f?Ni  =  Ry sz),  R1.5  (Rui  =  0.67Rysz),  R2  (Rn i  =  0.5Rysz)  and  R2.5 
(RNi  =  0.4Rysz)  anode  systems,  respectively.  The  decrease  in  total 
pore  interfacial  area  is  primarily  due  to  the  temporal  evolution  of 


Ni-pore  interfaces.  By  examining  the  evolution  of  pore  space,  it  can 
be  concluded  that  the  pores  are  growing  in  size  [Fig.  9(a)],  and 
evolving  into  a  more  equiaxed  shapes,  or  becoming  more  round, 
thus  minimizing  total  interfacial  area  [Fig.  9(b)].  Our  simulation 
results  for  the  decrease  in  total  pore  interfacial  area  are  in  good 
agreement  with  annealing  experiments  on  Ni-YSZ  porous  anodes 
with  Ni  and  YSZ  particles  of  similar  size,  where  a  —25%  decrease  in 
total  pore  interfacial  area  was  observed  [26]. 

Other  effective  materials  properties  are  also  extracted  during 
the  evolution  of  SOFC  anode  systems.  Fig.  10(a)  displays  a  plot  of  the 
pore  space  tortuosity  factor  a^,  which  affects  gas  transport 
mechanisms.  Initially,  anode  systems  with  smaller  initial  Ni  parti¬ 
cles  are  characterized  by  larger  values  for  the  tortuosity  factor.  A 
large  value  indicates  highly  twisted  paths,  which  can  slow  down 
fuel  gas  diffusion  processes.  Ionic/electronic  transport  mechanisms 
are  influenced  by  changes  in  effective  conductivities  during  the 
microstructural  evolution  of  SOFC  anodes.  Due  to  the  assumed 
static,  or  non-evolving,  nature  of  the  YSZ  phase,  the  effective  ionic 
conductivity  of  the  YSZ  phase  (7§fz  is  constant  and  approximately 
the  same  for  all  anode  systems.  Fig.  10(b)  in  turn  quantifies  the 
temporal  evolution  of  the  effective  electronic  conductivity  of  the  Ni 
phase  normalized  by  the  intrinsic  Ni  phase  bulk  conductivity 
g0.  Initially,  systems  with  smaller  Ni  particles  are  characterized  by 
larger  values,  but  as  these  systems  evolve  over  time,  they 
experience  the  largest  drop  in  <j^f. 

Finally,  we  examine  the  temporal  evolution  of  Ni  phase  perco¬ 
lation  degree  Pd.  Fig.  11(a)  is  a  plot  of  Ni  phase  percolation  degree  as 
a  function  of  coarsening  time.  It  can  be  clearly  seen,  that  during  Ni 
coarsening,  systems  with  smaller  initial  Ni  particles  experience  the 
largest  drop  in  Ni  phase  contiguity,  i.e.,  decrease  in  percolation 
degree.  A  close  examination  of  the  Ni  phase  grain  morphology  re¬ 
veals  that  the  loss  of  percolation  degree  may  be  associated  with  the 
presence  of  Ni  grain  boundaries.  The  two  panels  in  Fig.  11(b)  depict 
a  section  of  a  Ni  skeleton  in  the  R2.5  anode  system,  where  the  top 
panel  (early  time)  shows  a  connected  Ni  cluster.  In  the  bottom 
panel  (late  time)  of  Fig.  11(b),  it  can  be  seen  that  a  pinch-off  effect 
along  Ni  grain  boundaries  has  caused  the  cluster  to  lose  contiguity. 
In  an  experimental  study  by  Nelson  et  al.  [25]  on  SOFC  Ni-YSZ 
porous  anodes  under  operating  conditions  with  initial  Ni  and  YSZ 
particles  of  similar  size,  they  observed  trends  for  the  percolation 
degree  that  are  similar  to  the  ones  depicted  in  Fig.  11(a),  where  a 
continuous  decrease  in  Ni  phase  contiguity  was  observed. 

In  more  quantitative  terms,  the  temporal  evolution  of  Ni  phase 
percolation  degree  Pd  can  be  cast  in  a  functional  form  as  follows: 

PD(t)  -  PD(0)  =  Gre[(t/t*),  (13) 


Fig.  9.  Temporal  evolution  of  the  (a)  average  pore  space  radius  Rpore  normalized  by  the 
YSZ  particle  size  RY sz  and  (b)  total  pore  space  interfacial  area  (with  both  Ni  and  YSZ) 
Apore  normalized  by  the  initial  one  Apore(0).  Note  that  the  YSZ  phase  is  a  non-evolving 
one  and  its  particles  have  the  same  size  in  all  systems. 


-  R1  -  R1 


Fig.  10.  Temporal  evolution  of  (a)  pore  phase  tortuosity  factor  a*  and  (b)  Ni  phase 
electron  conductivity  normalized  by  the  intrinsic  Ni  bulk  conductivity  a0. 
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Fig.  11.  (a)  Temporal  evolution  of  Ni  phase  percolation  degree  PD,  which  is  a  measure 
of  percolating  Ni  clusters  that  span  the  simulation  box.  (b)  For  the  R2.5  anode  system, 
two  panels;  top  (early  time)  and  bottom  (late  time)  showing  a  section  of  a  Ni  skeleton. 
Circles  colored  in  red  show  a  pinch-off  effect  at  Ni  grain  boundary  regions  that  leads  to 
morphological  instability.  Grain  boundaries  are  colored  in  black.  (For  interpretation  of 
the  references  to  color  in  this  figure  legend,  the  reader  is  referred  to  the  web  version  of 
this  article.) 


where  t*  =  t  *  (Rm/Rysz)  denotes  a  characteristic  time  scale  that 
controls  the  loss  of  percolation  degree,  while  Gref  is  a  reference 
functional  form  for  the  percolation  degree  Pd.  The  self-similar 
behavior  of  Ni  phase  percolation  degree  is  depicted  in  Fig.  12(a), 
while  the  dependence  of  t*  on  the  ratio  Pni/^ysz  is  in  turn  depicted 
in  Fig.  12(b).  The  data  can  be  reasonably  well  fit  with  a  simple 
power  law  function  of  the  form  t*  —  (Pni/%sz)x*  where  %  —  4.1  ±  2.2 
denotes  the  effective  scaling  exponent. 

To  rationalize  this  strong  dependence  of  loss  of  percolation 
degree  on  initial  Ni  particle  size,  consider  the  following.  If  Ni  par¬ 
ticles  in  the  initial  stage  are  idealized  as  cylindrical  morphologies 
with  “beaded”  structures,  then  the  loss  of  Ni  phase  percolation 
degree  can  be  related  to  the  classic  Rayleigh  instability  of  inviscid 
liquid  jets  [51].  In  this  case,  the  cylindrical  liquid  jet  is  linearly 
unstable  to  perturbations  with  wavelengths  greater  than  the 
circumference  of  the  cylinder,  and  the  presence  of  such  perturba¬ 
tions  eventually  leads  to  the  breaking  of  the  jet  into  spherical 
droplets.  In  the  context  of  solid  cylinders,  McCallum  et  al.  [52] 
examined  analytically  the  stability  of  cylindrical  islands  on  rigid 
substrates  with  various  wetting  angles  and  demonstrated  that  the 
characteristic  time  for  the  onset  of  morphological  instability  for  the 


case  of  surface  diffusion  i s~L^/b,  where  b  is  a  parameter  that  en¬ 
compasses  various  materials  properties  such  as  diffusivity  and 
surface  tension,  and  Ld  is  a  length  scale  that  characterizes  the  cross- 
sectional  area  of  the  unperturbed  cylinder.  In  addition,  we  note  that 
grain  boundaries  also  play  an  important  role  in  the  loss  of  conti¬ 
guity  [cf.  Fig.  11(b)].  In  particular,  during  the  late  stages  of  the 
pinch-off  process,  grain  boundary  grooving  takes  place  [53].  By 
extending  Mullins’s  analysis  for  a  cylindrical  bicrystal,  it  can  be 
shown  that  the  time  to  pinch-off  displays  a  simple  R4  scaling  form, 
where  R  denotes  the  radius  of  the  unperturbed  cylinder.  Both  of 
these  simple  arguments  thus  suggest  a  t*  ~  ( Rm/Rysz )4  scaling  form 
for  the  characteristic  time  scale  for  loss  in  percolation  degree, 
consistent  with  our  numerical  simulations. 

4.2.  Electrochemical  performance 

In  order  to  examine  the  impact  of  Ni  coarsening  on  the  elec¬ 
trochemical  performance  of  SOFCs,  we  employ  the  microstructural 
features  and  effective  material  properties  that  previously  have  been 
extracted  from  evolving  anode  microstructures  as  input  parameters 
in  our  cell  level  model.  So  far,  no  absolute  physical  dimensions  of 
the  simulation  box  had  been  assigned  in  order  to  keep  the  results  as 
general  as  possible,  and  all  dimensions  were  scaled  with  respect  to 
the  characteristic  length  scale  d0.  For  the  electrochemical  analysis 
presented  in  this  section  we  assign  an  absolute  value  to  the  size  of 
the  non-evolving  YSZ  particles  and  choose  the  characteristic  length 
scale,  d0,  accordingly  in  order  to  obtain  physical  measures  of  the 
electrochemical  performance.  For  example,  a  YSZ  particle  size  of 
2  pm  (Rysz  =  1  pm)  corresponds  to  a  characteristic  length,  i.e.,  a 
voxel  size  of  d0  =  80  nm,  leading  to  a  (9.6  pm)3  sized  simulation 
box.  To  account  for  sample-to-sample  variations,  effective  material 
properties  are  averaged  over  four  ensembles  for  each  anode  system 
(Rl,  R1.5,  R2  and  R2.5). 

Power  density  vs.  current  density  plots  are  established  for  all 
SOFC  cells  during  the  microstructural  evolution  of  SOFC  anode 
systems.  Fig.  13(a)  shows  an  exemplary  plot  of  the  power  density  vs. 
current  density  for  a  SOFC  cell  having  an  R2.5  anode  system,  i.e., 
Rysz  =  1  pm  and  Rm/Rysz  =  0.4  at  the  initial  time  step,  t  =  0,  and 
after  a  coarsening  period  of  t  =  104t,  corresponding  to  ~360  h  in 
real  time  for  this  specific  configuration.  We  define  the  maximum 
power  density  as  the  peak  point  in  the  power  density  vs.  current 
density  curve.  A  decrease  of  the  maximum  achievable  power  den¬ 
sity  with  coarsening  time  can  be  observed  in  the  example,  indi¬ 
cating  that  the  changes  in  the  microstructure  caused  by  Ni 
coarsening  lead  to  a  performance  loss. 
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Fig.  12.  (a)  Self  similar  scaling  of  Ni  phase  percolation  degree  PD  for  all  SOFC  anode 
systems,  (b)  On  a  log-log  scale,  characteristic  time  f  =  t  *  (Fni/^ysz)  that  sets  the 
morphological  instability  vs.  Ni  to  YSZ  particle  size  ratio. 


Fig.  13.  For  a  YSZ  particle  size  of  2  pm  (RY sz  =  1  pm),  (a)  Power  density  vs.  current 
density  plots  for  a  SOFC  cell,  with  the  R2.5  anode  system  ( Rm/Rysz  =  0.4),  at  coarsening 
times  t  =  0  (initial)  and  t  =  104.  (b)  Maximum  power  density  vs.  coarsening  time  for  all 
considered  SOFC  anode  systems.  The  max.  power  density  is  defined  as  the  peak  point 
in  the  power  density  vs.  current  density  curve. 
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To  further  quantify  this  effect,  we  determined  the  temporal 
evolution  of  the  maximum  achievable  power  density  for  all 
considered  anode  configurations.  Taking  Rysz  =  1  pm,  the  role  of  the 
initial  anode  microstructure  on  the  evolution  of  the  electro¬ 
chemical  performance  during  Ni  coarsening  is  quantified  in 
Fig.  13(b).  In  particular,  initial  anode  microstructures  with  smaller 
Ni  particles  yield  higher  efficiencies,  i.e.,  values  for  the  maximum 
power  density,  which  can  be  explained  by  their  larger  amount  of 
available  reactive  TPB  [31].  This  is  in  agreement  with  recent 
experimental  findings  suggesting  that  SOFCs  with  anodes  having 
smaller  features  for  the  Ni  particles  yield  higher  maximum  power 
densities  than  the  ones  with  anodes  having  larger  Ni  particles  [54]. 

Due  to  the  microstructural  evolution  within  the  anode,  our  SOFC 
cell  model  predicts  a  drop  in  the  electrochemical  performance 
during  Ni  coarsening  for  microstructures  with  Rysz  =  1  pm.  This 
behavior  is  due  to  the  combined  effects  of  continuous  loss  of 
reactive  TPB  regions,  changes  of  the  electronic  conductivity  in  the 
Ni  phase  as  well  as  the  evolution  of  pore  sizes  and  pore  tortuosity 
factor.  While  anodes  with  initial  Ni  particles  smaller  than  YSZ  ones 
yield  a  higher  initial  maximum  power  density,  such  microstruc¬ 
tures  are  also  prone  to  the  strongest  decline  in  performance  as 
coarsening  takes  place.  The  drop  in  the  maximum  power  density 
becomes  smoother  as  the  initial  sizes  of  Ni  and  YSZ  particles 
become  comparable.  Note  that  SOFC  cells  with  R2.5  anode  systems 
(i.e.,  Rysz  =  1  pm  and  Rm/Rysz  =  0.4)  experience  a  sudden  drop  in 
the  maximum  power  density  at  late  times  t  >  7  x  103t.  This  is 
strongly  correlated  with  the  loss  of  reactive  TPB  regions  observed  in 
Fig.  7(b),  which  in  turn  is  associated  with  the  loss  of  the  percolation 
degree  Pd  in  the  Ni  phase  as  illustrated  in  Fig.  11(a). 

Next,  we  repeat  this  approach  for  a  variety  of  YSZ  particle  sizes 
in  a  range  from  Rysz  =  50  nm  to  Rysz  =  4  pm.  Thus,  we  determine  a 
corresponding  characteristic  length  scale  d0  for  each  specific  Rysz 
value.  The  purpose  of  this  study  is  to  examine  the  impact  of  the 
absolute  size  of  the  YSZ  particles  in  the  non-evolving  backbone 
structure  on  the  electrochemical  performance  during  coarsening. 
As  before,  for  each  value  of  Rysz,  we  consider  four  cases  with 
different  size  ratios  between  Ni  and  YSZ  particles.  Each  anode 
microstructure,  characterized  by  its  combination  of  Ni  and  YSZ 
particle  size  and  ratio,  shows  a  different  initial  maximum  power 
density  before  coarsening  takes  place.  This  is  caused  by  the 
different  scaling  behavior  of  the  effective  material  parameters 
used  as  input  parameters  for  our  cell  level  model:  for  example,  at  a 
given  volume  fraction  for  the  anode  phases,  scaling  a  given 
morphology  to  half  of  its  original  size,  i.e.,  reducing  the  size  of  the 
YSZ  particles  (and  also  the  Ni  particles)  by  50%,  consequently  leads 
to  a  reduction  of  the  pore  size  of  50%,  while  the  TPB  density  is 
increased  by  a  factor  of  four  [31].  Since  we  are  interested  in 
monitoring  relative  changes  of  the  performance  caused  by  the 
microstructural  evolution  during  Ni  coarsening,  we  normalize  the 
values  found  for  the  maximum  power  density  in  coarsened  con¬ 
figurations  by  their  respective  initial  ones.  Contour  maps  of  the 
normalized  power  density  are  generated  in  order  to  visualize  the 
results  and  identify  regions  of  initial  YSZ  and  Ni  particle  sizes  for 
which  the  performance  is  potentially  stable  against  coarsening  of 
the  Ni  phase. 

Fig.  14(a)  and  (b)  shows  such  maps  of  the  normalized  maximum 
power  density  at  (non-dimensional)  simulation  times  of  1.0  x  103 
and  1.0  x  104,  respectively.  [Note  that  different  values  of  Rysz  are 
associated  with  different  values  of  the  characteristic  length  scale  d0 
as  mentioned  before.  According  to  Eq.  (5),  this  leads  in  turn  to 
different  values  of  the  characteristic  time  t  as  Rysz  is  varied.]  These 
maps  reveal  that  the  size  of  the  YSZ  particles  in  the  non-evolving 
backbone  structure  as  well  as  the  relative  size  of  the  Ni  particles 
have  a  strong  influence  on  how  the  performance  is  impacted  by  Ni 
coarsening,  as  will  be  discussed  next. 
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Fig.  14.  Maps  of  the  maximum  power  density  normalized  by  its  initial  value  for 
various  combinations  of  YSZ  particle  size  and  Ni  to  YSZ  particle  size  ratio  ( RmIRysz )  at 
coarsening  times  (a)  t  =  1.0  x  103  (b)  t  =  1.0  x  104. 

A  close  examination  of  the  maps  in  Fig.  14(a)  and  (b)  reveals  that 
coarsening  of  the  Ni  phase  causes  a  degradation  of  electrochemical 
performance  for  YSZ  particles  of  size  larger  than  approximately 
Rysz  >  0.75  pm.  Furthermore,  for  systems  with  smaller  YSZ  particles 
(Rysz  ^  0.5  pm  and  %i  <  Rysz)  and  for  the  simulation  times  attained 
in  this  study,  there  exists  a  time  window  within  which  Ni  coars¬ 
ening  can  even  enhance  the  electrochemical  performance.  When 
comparing  both  maps,  i.e.,  maps  of  the  relative  performance  change 
obtained  at  different  simulation  time  steps,  both  observed  trends 
become  more  pronounced  at  longer  simulation  times. 

We  suggest  the  following  explanation  for  these  results.  First,  we 
note  that  porous  anode  microstructures  with  a  relatively  large  YSZ 
particle  size,  i.e.,  Rysz  >  0-75  pm,  also  have  a  sufficiently  large 
average  pore  size,  so  that  the  gas  transport  through  the  pore 
network  to  and  from  TPB  sites  can  easily  be  sustained.  Flowever,  the 
amount  of  active  TPB  represents  the  limiting  factor  for  the  elec¬ 
trochemical  performance  of  such  electrodes.  Ni  coarsening  leads  to 
a  reduction  of  the  amount  of  TPB  sites  (cf.  Fig.  7),  so  that  the  overall 
electrochemical  performance  degrades  as  Ni  coarsening  takes 
place. 

In  comparison,  anode  microstructures  with  a  relatively  small 
YSZ  particle  size,  i.e.,  Rysz  ^  0.5  pm,  and  Ni  particles  of  even  smaller 
size,  exhibit  a  sufficient  amount  of  TPB,  so  that  the  number  of 
reactive  sites  does  not  become  rate  limiting.  Flowever,  the  small 
diameters  of  the  porous  network  become  more  obstructive  for  the 
diffusive  gas  transport,  especially  for  pore  sizes  where  the  Knudsen 
diffusion  becomes  comparable  to  or  even  dominant  over  the  mo¬ 
lecular  gas  diffusion.  We  capture  this  effect  in  our  SOFC  cell  level 
model  by  using  the  dusty-gas  approach  that  takes  both  diffusion 
mechanisms  into  account.  When  Ni  coarsening  takes  place  in  such 
electrodes,  the  gas  channels  in  the  pore  network  widen  and  open 
up  [cf.  Fig.  9(a)],  which  yields  improved  and  alleviated  conditions 
for  the  diffusive  gas  transport.  While  Ni  coarsening  again  leads  to  a 
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reduction  of  TPB,  this  effect  now  becomes  balanced  or  even  over¬ 
compensated  by  improved  gas  transport  properties  due  to 
widening  pore  sizes.  Note  that  the  electronic  conductivity  of  the  Ni 
phase  also  evolves  during  coarsening  [cf.  Fig.  10(b)].  However,  as 
long  as  a  reasonable  extent  of  Ni  phase  contiguity  is  maintained,  we 
judge  this  to  have  only  a  minor  influence  on  the  performance,  since 
the  ionic  conductivity  of  the  YSZ  phase  is  orders  of  magnitude 
lower  than  the  electronic  conductivity  of  Ni. 

In  terms  of  electrochemical  efficiency,  our  model  predicts  that 
anodes  in  the  latter  regime,  i.e.,  for  Rysz  <  0.5  pm  and  R^i  <  Rysz  are 
better  capable  of  withstanding  the  microstructural  changes  asso¬ 
ciated  with  Ni  coarsening,  provided  that  the  YSZ  backbone  does  not 
evolve  with  time  and  the  Ni  phase  maintains  a  reasonable  level  of 
percolation  degree,  or  contiguity.  Therefore,  we  propose  that  tun¬ 
ing  the  initial  particle  sizes  of  Ni  and  YSZ  in  anode  cermets  can  have 
a  crucial  impact  on  the  coarsening  behavior  and  the  closely  linked 
temporal  evolution  of  the  electrochemical  performance  during  fuel 
cell  operation. 

5.  Conclusions 

In  summary,  an  integrated  modeling  approach  was  presented  in 
order  to  quantify  the  effect  of  microstructural  evolution  of  Ni-YSZ 
cermet  anodes  on  the  electrochemical  performance  of  SOFCs.  A 
mesoscale  phase-field  model  was  developed  that  captures  Ni 
coarsening  in  Ni-YSZ  cermet  anodes,  accounts  for  the  poly¬ 
crystalline  nature  of  Ni,  and  tracks  the  evolution  of  several  micro- 
structural  features  that  influence  the  electrochemical  performance. 
These  features  were  then  extracted  and  used  as  inputs  in  a 
macroscale  cell  level  model  that  encompasses  gas  transport,  elec¬ 
tron  and  ion  conduction  processes  as  well  as  the  electrochemical 
reactions  in  order  to  examine  the  role  of  Ni  coarsening  on  elec¬ 
trochemical  performance.  While  the  electrochemical  model  itself  is 
time-independent,  its  inputs,  which  represent  effective  properties, 
were  extracted  from  the  phase  field  model  as  Ni-YSZ  porous  an¬ 
odes  evolved  over  time.  Therefore,  the  cell  level  model  is  capable  of 
capturing  the  temporal  evolution  of  the  electrochemical  perfor¬ 
mance  of  SOFCs.  Initial  anode  microstructures  were  idealized  as 
binary  mixtures  of  spherical  Ni  and  YSZ  particles  and  the  integrated 
modeling  approach  was  employed,  in  a  parametric  study,  to 
examine  the  role  of  initial  anode  microstructure  on  the  morpho¬ 
logical  stability  of  cermet  anodes  and  electrochemical  performance 
of  SOFCs. 

Our  simulation  results  revealed  simple  scaling  laws  for  the 
decrease  in  overall  TPB  length  and  loss  in  Ni  percolation  degree,  or 
contiguity.  A  Ni  phase  size  effect  was  revealed,  such  that  systems 
with  initially  smaller  Ni  particles  maintained  higher  density  of  TPBs 
but  lost  contiguity  faster  than  the  ones  with  larger  Ni  particles. 
Finally,  the  overall  electrochemical  performance  displayed  a  non¬ 
monotonic  time-dependence  during  coarsening. 

For  initially  coarse  microstructures,  the  interplay  between  the 
ceaseless  decrease  in  reactive  TPB  density  and  the  increase  in 
average  pore  radius  led  to  degradation  in  electrochemical  effi¬ 
ciency.  On  the  other  hand,  for  initially  fine  microstructures,  these 
competing  factors,  i.e.,  the  decrease  in  TPB  density  and  rapid  in¬ 
crease  in  average  pore  space  radius  led  to  stabilizing  the  degrada¬ 
tion  of  electrochemical  efficiency.  This  could  be  an  indication  that 
at  these  length  scales,  electrochemical  performance  is  dominated 
by  gaseous  transport  processes  through  the  pore  phase. 

In  broader  terms,  our  approach  highlights  the  importance  of  the 
initial  anode  microstructure  on  both  the  morphological  stability  of 
Ni  phase  and  electrochemical  performance  of  SOFCs.  Several 
physical  phenomena  are  identified  as  possible  mechanisms  that 
control  Ni  phase  stability  and  coarsening  rates  in  SOFC  porous 
anodes.  Our  integrated  modeling  approach,  which  is  informed  by 


the  local  microstructure,  can  be  used  as  a  design  tool  to:  i)  Identify 
combinations  of  morphological  parameters  that  yield  anode  sys¬ 
tems  with  optimal  characteristics  for  the  microstructural  stability 
and  electrochemical  performance  of  SOFCs,  ii)  utilize  the  trends  and 
scaling  laws  for  TPB,  Ni  contiguity  and  pore  radius  to  predict  the 
performance  of  SOFCs  under  extended  use,  and  iii)  explore  other 
avenues  aiming  at  stabilizing  Ni  phase  coarsening,  such  as  doping 
Ni  to  increase  the  wettability  with  the  YSZ  backbone. 

Finally,  the  role  of  electrochemical  reactions  on  the  evolution  of 
several  materials  properties  that  affect  Ni  phase  coarsening  dy¬ 
namics  is  not  accounted  for  in  our  approach.  For  example,  Ni  phase 
atomic  mobility  depends  on  the  local  temperature,  which  is  influ¬ 
enced  by  electrochemical  reactions.  Further  improvement  to  our 
modeling  approach,  which  is  the  subject  of  a  future  research 
exploration,  would  incorporate  the  effects  of  electrochemistry  on 
Ni  phase  coarsening  dynamics  in  SOFC  anodes.  In  addition,  the  role 
of  coarsening  on  the  redox  stability,  alternating  oxidizing  and 
reducing  conditions,  is  an  important  factor  in  the  design  of  durable 
SOFC  anode  microstructures.  To  this  end,  we  utilized  a  recently 
developed  model  [55]  in  combination  with  the  coarsening 
approach  presented  in  this  work  to  examine  the  role  of  Ni  phase 
microstructural  evolution  on  the  redox  stability  of  SOFC  cermet 
anodes  also  consisting  of  binary  mixtures  of  initially  spherical  Ni 
and  YSZ  particles,  the  results  of  which  are  published  elsewhere 
[62]. 
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